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Abstract 

We use the hyperbolic subdiflFusion equation with fractional time derivatives (the generalized 

j^ , Cattaneo equation) to study the transport process of electrolytes in media where subdiffusion 

occurs. In this model the flux is delayed in a non-zero time with respect to the concentration 

gradient. In particular, we obtain the formula of electrochemical subdiffusive impedance of a 

^ . spatially limited sample in the limit of large and of small pulsation of the electric field. The 

boundary condition at the external wall of the sample are taken in the general form as a linear 

combination of subdiffusive flux and concentration of the transported particles. We also discuss 

'■^ ■ the influence of the equation parameters (the subdiffusion parameter and the delay time) on the 

^ ' Nyquist impedance plots. 

O' 

o 

1 Introduction 

^ ' Subdiffusion occurs in systems where mobility of particles is significantly hindered due to internal 

CO , structure of the medium, as in porous media, gels or amorphous semiconductors [HE]- The subdiffu- 

^^ ' sion is characterized by a time dependence of the mean square displacement of transported particle 

r^ , (Ax^) = YTi+S^y where Da is the subdiffusion coefficient measured in the units m^/s" and a is the 

• ' subdiffusion parameter of a value within < a < 1 range. For a = 1 one deals with the normal 

J— ^ . diffusion. 

QQ ' The subdiffusion has been recently extensively studied. While the phenomenon is theoretically 

(^ . rather well understood, there are very few reported experimental investigations (e.g. [2])- The 

method of impedance spectroscopy was used to experimentally study subdiffusion in porous media 
such as nanopore electrode [$1, cement [11 [H [6], tooth enamel [7] and gels [8]. The theoretical 
analysis of subdiffusion impedance was presented by f9] who used the following parabolic subdiffusion 
equation (PSE) with fractional time derivative. The equation reads as 



X 



dt "ati-" 9a;2 ' ^^> 

where the Riemann-Liouville fractional time derivative, which is defined for a > as [TT| [TO] 

d"f{t) _ 1 5" /■*^^,_ f{t') 



dt°' T{n-a)dt" Jo {t-t'Y+°'-"- 

For a — 1 Eq. ^ converts into the normal diffusion equation. 

For the initial condition C(x,0) — 5{x), where 5 is the Dirac-delta function, the solution of 
Eq. ([1]) (the Green's function) has non-zero values for any x and t (i > 0). Thus, even for small 



times, a finite amount of the substance exists at very large distances from the origin, what can 
be interpreted as infinite propagation velocity of some of the diffusing particles. To avoid this 
'unphysical property' Cattaneo derived the hyperbolic normal diffusion equation which Green's 
function has non-zero values for finite x |12[ [T3] . The equation is based on the assumption, that 
the flux is delayd by time period r with respect to the concentration gradient. For many 'typical 
systems' (as the membrane one) it is hard to observe the difference between the solutions of parabolic 
and hyperbolic (sub)diffusion equations even for relatively large values of r p3]. However, in some 
processes the non-zero parameter r plays crucial role. The example is the diffusion in a system 
where boundary conditions are given by functions quickly changing in time. Such a situation occurs 
in the electrochemical system with (sub)diffusion impedance. As far as we know, the Cattaneo 
equation was used to study electrochemical impedance only for a system where normal diffusion 
occurs [15, 16j, except our work [17J where the subdiffusion impedance was considered in a system 
with fully absorbing wall. Articles pubHshed so far mostly concentrated on homogeneous systems, 
however, a more complex system containing few diffusion layers was also studied in [18| . 

In this paper we present a theoretical foundation for studies of subdiffusion impedance using a 
hyperbolic equation. We apply the hyperbolic Cattaneo equation with the fractional time derivative 
to model the subdiffusion impedance of a homogeneous sample of finite thickness, where the bound- 
ary condition at the sample surface is assumed as Hnear combination of fiux and concentration. We 
find an infiuence of the parameters a and r on the final formula describing the impedance of the 
subdiffusive medium, particularly for high and for low ac-voltage frequency. 

2 The generalized Cattaneo equation 

The phenomenological derivation of Cattaneo equation is based on the assumption that the fiux of 
the particles J is not generated by the concentration gradient instantaneously (as in the process 
described by parabolic diffusion equation), but it is delayed in time by r [T2j, what is refiected by 
the following relation 

J(x,t + T) = -I?„^fe^. (2) 

Expanding the left-hand side of Eq. ^ into power series with respect to r and assuming that the 
parameter r is sufficiently small, one gets 

Combining Eq. |[3)) with the continuity equation 

dC{x,t) dJ{x,t) 
dt dx 

one obtains the hyperbolic normal diffusion equation 

d^ C{x,t) dC{x,t) _ d^C{x,t) 

dt^ ^ dt " " dx^ • ^' 



(4) 
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The parabolic subdiffusion equation can be derived form Continuous Time Random Walk formaHsm 
or using the phenomenological approach. In the letter case one sets the Riemann-Liouville fractional 
time derivative of the order 1 — a in the right-hand side of Eq. ^ (with r = 0) or replaces the time 
derivative of the first order in Eq. Q to Caputo fractional time derivative of the order a. In similar 
way one can obtain the hyperbolic subdiffusion equation. In the paper [13] there were proposed three 
possible generalizations of the Cattaneo equation ^ to the one with the fractional time derivative, 
each one supported by a different scheme. In each of them the fractional time derivative replaces 
the derivative of natural order in Eqs ([3]) or |[4l), or it is added to the right-hand side of Eq. ^. 



The schemes provide a different forms of hyperboUc subdiffusion equation, which are not equivalent 
of each other. In the next section we use the last scheme - which is the most natural in our opinion 
- to derive the hyperbolic subdiffusion equation. Thus, the generalized Fick equation is given by 



J{x,t + T) 



-Dr, 



91-" dC{x,t) 



'dt^ 



dx 



(6) 



Equation ^ ensures that changes of the flux due to the concentration gradient are delayed by the 
time r. Assuming that r <C i and keeping linear terms with respect to r in the series expansion of 
l.h.s. of Eq. (ini), we get 



j(.,t)+.^4:^^-i..-^""^^(-'*) 



dt ""ati-" dx ■ 

Combining Eq. Q with the continuity equation Q, we obtain the generalized Cattaneo equation 



(7) 



d'^C{x,t) dC{x,t) 91-" 92(^(x,t) 



dt^ 



dt 



dt^ 



dx^ 



(8) 



3 The system 



Let us assume that at a; = there is the oscillating overvoltage ri{t) = Esm{ujt) which causes the 
oscillation of the concentration on the surface layer according to the formula 



where eq denotes a derivative computed in the local equilibrium. Thus, we have 

C(0,t) =Cosin(tji), 



(9) 



(10) 



where Co = RwQ^E, Rw = ^ ( dc ) • '^^^ conduction current I{t) at a; = corresponds to the 

flux of diffusing particles J(0, t) 

I{t) ^ qAJ{0,t), (11) 

where q is the charge of diffusing particle and A is the area of the sample surface. The surface located 
at a; = L can be treated as fully absorbing, partially absorbing or fully reflecting wall (FigH]). 



r\(t)=EffSin((i)t) 
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Figure 1: The system under considerations, 77 denotes the overvoltage, E - its amplitude. 

The second boundary condition, fixed at the wall x — L, can be chosen in different ways de- 
pending on the properties of the wall. Usually, for the system described by parabolic (sub)diffusion 
equation one can adopt the boundary conditions as follows: 



For the fully absorbing wall one gets [19] 

C{L,t)^0. (12) 

For the fully reflecting wall there is [19] 

J(L,i) = 0. (13) 

For the partially absorbing wall the particle absorbed by the wall cannot return to the system 
and the boundary condition is given by |20[ [21] 



,,, , dC(x,t) 



dx 



(14) 



The boundary condition for the fully absorbing or fully reflecting wall was studied in majority works 
[22l E], however, the radiation boundary condition lfT4| was also used [HI O El]- In this work 
we assume that the boundary condition at a; = L is given in a general form and it is a linear 
combination of the flux and concentration, where the flux is delayed in time by r with respect to 
the concentration. So, the boundary condition is 

aLJ{L,t + T)+bLC{L,t)^Q. (15) 

4 Diffusion impedance 

The impedance of electrochemical system Z{iuo) can be defined as its response to a voltage or current 
perturbation from a steady-state situation |22[ [25] 

T){iuj) 

^M = 77--' (16) 

where f}{iLj) and I{iuj) are the Laplace transforms of the overvoltage and current perturbation, uj is 
the angular frequency uj = 27r/, where / is frequency. A plot of the real and imaginary parts of the 
impedance in the complex plane {ReZ, —ImZ) as the frequency (treated as a parameter) is swept 
over a given range is called the Nyquist plot. 

The impedance of the diffusion layer is called the Warburg impedance. For the layer of the 
infinite thickness the impedance is defined as 

Zi^u;)^^^^il~^), (17) 



where R is the diffusion resistance. On the Nyquist plot the Warburg impedance is presented by the 
straight half line with the slop angle 7r/4 passing through the origin of coordinates. In real systems 
the diffusion layer has a finite thickness. Let us assume, that the diffusion layer is bordered by planes 
localized at x = and x = L. The perturbation of the voltage is applied to the medium at a; = 0. 
The characteristic angular frequency is defined as 



D 

Z2 



^d = T2 , (18) 



where D is the diffusion coefficient. The frequency (fTSl) is proportional to the inverse of the average 
time necessary for an ion to cross the sample thickness. When w ^ w^ the size of the sample plays 
no role in the ion diffusion and the impedance is the Warburg impedance (fl7|) . However, for a 



low frequency the ions can be absorbed by the opposite wall before they change direction of their 
movement. 

From Eqs. (fT6|) . |[9]) and (fTTI) one obtains the relation 



Z{ioj) = R 



w 



J{0,iuj) 



(19) 



5 SubdifFusion impedance 

We assume that the transport process is described by the Cattaneo equation ^ with the following 
initial conditions 



dt 



0. 



t=o 



The Laplace transform of ^ for the initial conditions l(20|) is 

Ts^C{x, s) + sC{x, s) - i?„si-"^^^^^. 

The general solution of Eq. (|2T|) is 

C{x, s) — Bi exp [7(s)a;] + B2 exp [—j{s)x] . 

where 

lis) 

The Laplace transform of the flux is 






--VT 



J{x, s) — —Da 



s^^" dC{x,s) 



1 + Ts dx 



Combinig Eqs. (dH), ^ and (I22l)-(l2ll) we obtain 



Z{s) = Rw 



m 



bL sinh(7(s)L) — aLX[s) cosh(7(s)L) 



6l cosh(7(s)i) — aiA(s) sinh(7(s)L) 



where 



A(.) = .i-"/^y^ 



Dc 



For the layer with infinite thickness {L — > 00) the impedance (|25l) has the following form 



Z{iijj) ~ Rw 



\/l + TIOJ 



(20) 

(21) 

(22) 
(23) 

(24) 

(25) 
(26) 

(27) 



For T = and a = 1 Eq. (|27|) corresponds to the classical Warburg impedance [22]. For subdiffusive 
systems the relation on the impedance (fTSJ) should be repleaced by ujd = {D/L'^y/°' [9]. 

When w ^ cx), a substantial infiuence of r can be inferred form Eqs. (|25| . When w ^ 00 from 
Eqs. ^ and ([26]) we obtain 



• For T 7^ 



Z{iuj) 






1 -a 



1-a 



(28) 



and the Nyquist plot of the impedance is a linear function passing through the origin of 
coordinates with the angle slope ip given by 



1 



Lp = -K- 



Let us note that for < a < 1 there is (^ G (0, 7r/4). 



• For T = 



where 



Z(iijj) 



Rw \ ( l-a/2\ . . / l-a/2 

l-a/2 
Lp = n 



(29) 

(30) 
(31) 



and ip e (7r/4, 7r/2). 
For low uj one obtains 



Z{iLu) = — 



Rw 

Xiito) 



aiA(icj) + bLLj{iuj) — aLL\{iuj)^{iuj) 



(32) 



An analysis of Eq. l(32|) leads to the following conclusions. For cj — > oo and a e (0, 1) that the 
slope of the plot is (/? = 7r(l — a)/2 when Bl ^ (for the partially or fully absorbing wall) and 
ip = 7r/2 when Bl — (for reflecting wall). We note that for oj <^ l/r the terms which contain the 
T parameter can be neglected in the above formulas. Thus, in all cases (p is also independent of t. 

Calculating ReZ and ImZ from l|25p for s — iuj we obtain the Nyquist plots (Figs. [2H5l) with 
several values of the parameters r and a. Our calculations were done for uj G (10~^, 10^), Rw = 1, 
L = I and Da = 1 (all quantities are in arbitrary units) . For larger values of uj the points on the 
plots are located near the origin and for larger r the curves are located near the ReZ axis. 
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Figure 2: The Nyquist plots for a = 0.4 and for r given in the legend. 

On Figs. ^ and ^ the plots for t = and t — 0.01 are practically undistinguishable, but for 
larger values of the subdiffusion parameter a these plots differ from each other. The curves with 
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Figure 3: The Nyquist plots for a. = 0.6. 
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Figure 4: The Nyquist plots for a = 0.8. 



larger values of the t parameter lie nearer the axis ReZ . The plots suggest that the curves for 
different values of r converge into one curve for very small values of u. For relatively large values of 
T the Nyquist plots show 'chaotic' behavior, which is stronger when a increases. Therefore, in the 
presented cases we did not consider values of r larger than the one on plots in Figs. (121) - Q. 
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Figure 5: The Nyquist plots for the normal diffusion case (a = 1.0). 



6 Final remarks 

The main result of our work is Eq. l(25l) with the asymptotic formulas (|28l) - l(32|) . These functions 
illustrated by the plots Figs. EHU show the following: 

For Lj ^ oo the Nyquist plot is the linear function passing through the origin of coordinates 
with the slop angle Lp dependent on subdiflusion parameter a and given by Eq. (|29)) and (pT]) . It 
is interesting that Lp does not depend on the parameter r explicitly, but the dependence is 'hidden' 
and is manifested by the interval to which the parameter belongs. For a given a G (0, 1) and t ^Q 
if> = 7r(l — q;)/4 < 7r/4 and ip = 7r(l — a/2)/2 > 7r/4 when r = 0. For non-zero r, ip does not 
dependent on r and for all values of (p it is independent of the boundary conditions at x — L. In the 
paper [24j the case oi ip y^ tt/A was interpreted somewhat different, mainly as a presence of so-called 
'constant phase element' (CPE) in the system. 

For relatively high w and for t ^^ the Nyquist plots show the 'chaotic' and 'oscillating' behavior, 
increases with increase of r and a parameters. The physical interpretation of this fact can be as 
follows. Namely, the periodic changes of concentration at the surface a; = generated the flux 
which is delayed in time by r with respect to the concentration gradient. When the oscillations of 
concentration are very rapid, the flux does not keep up with concentration changes, so the additional 
factor contributing to the total impedance is created. The difficulties in movement of the ions increase 
when a decreases causing a decrease in the flux. Thus, the r parameter less influences the transport 
process when a is smaller. Such a behavior is observed on the presented plots, where the curves 
obtained for t — and for t = 0.1 differ slightly form each other when a = 0.4 and a — 0.6, but 
the difference is relatively large for a — 0.8 and a ~ 1. For low lo the plots are dependent on the 
boundary conditions. According to the Eq. (|32| . for fully reflecting wall and u ^ l/r one gets 
ip — 7r/2, which is independent of a. For fully or partially permeable wall we observe that the plots 
become the linear where ip = 7r(l — a)/2. 

As we mentioned in the Introduction, there are a few methods to extract the value of subdiffusion 
parameters from experimental data. The considerations presented in this paper show that it is 
possible to determine the values of parameters of the system from the Nyquist plots obtained from 



experimental data. Measuring the slope angle for high w, one can extract the subdiffusion parameter 
a according to the formulas ((29l) and l(3T|) . Such a method was used to extract the parameter a 
for the lithium transport in gel electrode [8j. The authors found that Lp = 32°, it interpreted the 
results as subdiflusive transport of the particles inside the electrode. Since ip < 45°, we deduce that 
the transport studied in [8| can be described by the Cattaneo equation with non-zero parameter r. 
Unfortunately, the value of slope angle is not sufficient to extract r from experimental data. To find 
this parameter one should perform more detailed studies where the parameter representation of the 
Nyquist plots is taken into considerations. 

For the system described by parabolic subdiffusion equation ^ one obtains ip > tt/A. To achieve 
the situation where ip < 7r/4, the subdiffusion equation different from ^ was involved into the model 
[9] . Namely, it was assumed that the subdiffusion is described by one of the following equations 

dt^-^ " " dx^ ' ^ > 

when both of the equations contain the Riemann-Liouville fractional derivative. However, physical 
meaning of Eqs. ((33l) and (|34l) is rather unknown since they were derived only on the phenomeno- 
logical way, where the time derivative of natural order was replaced by the fractional one in the 
continuity equation and/or in the Pick's law. These equations were not derived on the base of 'mi- 
croscopic' model such as Continuous Time Random Walk formalism. We note that the Eq. J!]) is 
equivalent to l(33|) if in the latter one the Riemann-Liouville fractional derivative is replaced by the 
Caputo derivative. In our paper we show that the slope ip < 7r/4 is achieved form the model based 
on hyperboHc subdiffusion equation, which has clear physical interpretation. 
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